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ABSTRACT 

We present the observations of the starburst galaxy M82 taken with the Herschel SPIRE Fourier Transform Spectrometer. The spectrum (194-67 1 
yum) shows a prominent CO rotational ladder from J = 4-3 to 13-12 emitted by the central region of M82. The fundamental properties of the gas 
are well constrained by the high / lines observed for the first time. Radiative transfer modeling of these high-S/N 12 CO and 13 CO lines strongly 
indicates a very warm molecular gas component at ~ 500 K and pressure of ~ 3 x 10 6 K cm" 3 , in good agreement with the IT rotational lines 
measurements from Spitzer and ISO. We suggest that this warm gas is heated by dissipation of turbulence in the interstellar medium (ISM) rather 
than X-rays or UV flux from the straburst. This paper illustrates the promise of the SPIRE FTS for the study of the ISM of nearby galaxies. 
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1. Introduction 

Starburst galaxies provide us with the opportunity to study star 
formation and its effect on the interstellar medium (ISM) in ex- 
treme environments. These galaxies combine large central gas 
concentrations and high ionizing radiation fields, resulting in 
bright molecular, neutral and ionized gas emission lines. 

At a distance of 3.9 Mpc (ISakai & Madordl999h . M82 is the 
most well-studied starburst galaxy in the local universe, and it is 
widely used as a starburst prototype in cosmological stud ies. Its 
infrared luminosity (5.6 x 10 10 L^. ISanders et al.l|2003l) corre- 
sponds to a star-formation rate of 9.8 M Q yr _1 , which has almost 
certai nly been enhanc ed by its interaction with M8 1 and NGC 
3077 (lYun et aDll993h. With a re ported molecular gas content 
of 1.3 x 10 9 M Q dWalter et all2002h . its bright emission lines of 
CO and other molecules allow us to study its ISM in g reat detail 
(IShen& Loll 1991 IWalter et alJl200l IWard et al.ll200l . 

Far-infrared fine structure lines were used to constrain the 
physical properties of t he ionized gas and p hoto-dissociation re- 
gions (PDRs) in M82. IColbert et all d!999t) found that the ion- 
ized gas emission can be reproduced with a 3-5 Myr old instan- 
taneous starburst and a gas density of 250 cm 3 , while the PDR 
component is best fit with a density of 2000 cm 4 , in pressure 
equilibrium with the ionized phase. 



Stellar evolution an d photoionization models 
(Forster Schreiber et alj 120031) indicate a series of a few, 
Myr-duration starbursts with a peak of activity 10 Myr ago 
in the central regions, and 5 Myr ago in the circumnuclear 
ring. Models of the PDR and molecular emission as a set of 
non-interacting hot bubbles driving spherical shells of swept-up 
gas into a surrounding uniform medium also predict a starburst 
age of 5-1 Myr, bu t fail to match the observed far-infrared 
luminosity (Yao 2009]). 



* Herschel is an ESA space observatory with science instruments 
provided by European-led Principal Investigator consortia and with im- 
portant participation from NASA. 



The strengths of the CO lines place fundamen tal constraints 
on the physical properties of the molecular gas. iTilanus et al.l 
(1199 ll) fitted 12 CO and 13 CO lines from the central starburst up 
to / = 3-2 with a single-component model with temperatures 
of 30-55 K and densities of 3-7 xlO 3 cm 3 . IWild et all (Il992l) 
used lines up to the CO J = 6-5 transition to refine these pa- 
rameters to 40-50 K and ~ 10 4 cm 3 , while HCN and HCO + 
lines suggested densities greater than 3 x 10 5 cirT 3 are present. 
iPetitpas & Wilsonl (2000) showed evidence for a temperature or 
density gradient across the starburst region. I Weiss et all (120051) 
showed that CO emission up to J — 3-2 is dominated by more 
extended regions while higher / transitions originate in the cen- 
tral disk. 

In this paper, we present observations of M82 with Herschel 
dPilbratt et all 1201 Oh using the SPIR E Fourier Transform 
Spectrometer (FTS) (Griffin et al.1 120101) . which measures the 
complete far-infrared spectrum from 194 to 671 //m. This spec- 
tral region is particularly interesting for probing the peak of the 
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Fig. 1. Apodized spectrum of M82 corresponding to a 43'.'4 
beam, where red and blue lines represent data from the long- 
and the short-wavelength FTS bands respectively. Filled circles 
show SPIRE fluxes measured in the same beam. 



CO spectral line energy distribution (SLED) in gas-rich galax- 
ies. The wealth of lines across a continuous spectral region al- 
lows for unprecedented precision in modeling the physical and 
chemical properties of the molecular ISM. Here, we focus on the 
measurement and analysis of the CO rotational transitions from 
the central starburst in M82. 



2. Observations and data reduction 

The galaxy M82 was observed by the SPIRE FTS in the high 
spectral resolution (FWHM=0.048 cm -1 ), point-source mode, 
on 2009 September 21 as a performance verification target. The 
total integration time was 1 332 secon ds. The data w e re pro - 
cessed and calibrated as described in ISwinvard et a l. (2010). 
Only data from the central detectors in the two FTS bands are 
presented here. 

The beam size of the FTS bolomet ers varies with w avelength 
across the individual bands (see Swinyard et al. 2010), and the 
spatial extent of the M82 central starburst is comparable to the 
beam size (mean FWHM ~19" and 35" for the short- and the 
long-wavelength bands respectively). For a proper comparison 
with models, the spectrum must be scaled appropriately to a sin- 
gle beam size by a source-beam coupling factor (rj c (v)). This 
factor was obtained by convolving the M82 SPIRE photometer 
map at 250 yum (Roussel et al. 2010), which has a beam FWHM 
of 18'.'1 (Griffi n et alj|20101) . with appropriate Gaussian profiles 
to reproduce the light distribution as seen by FTS bolometers 
at different beam sizes. The value of rj c (v) is then given by the 
ratio of the beam-integrated flux density of the map convolved 
to the beam size corresponding to the given frequency (v) to the 
beam-integrated flux density of the map with the largest beam 
size (43'.'4); its values goes from 1 to 0.42. This implicitly as- 
sumes that the dust and CO emission distributions within the 
beam are the same at all frequencies. 

We opted to use the extended-source calibratecQ spectrum 
because the point-source calibration was more noisy and suf- 
fered from significant uncertainties below 600 GHz. We found, 
however, that the extended-source calibrated spectrum corrected 
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Fig. 2. 12 CO SLED for M82 for a 43'.'4 beam as meas ured in this 
work (open squares) with \a statistical error bars. IWard et all 
(2003) ground based data are shown by open circles. 



for source-beam coupling is around a factor of 2 fainter than 
photometry for the same beam. We thus scaled the spectrum 
to match the photometry in the three bands by applying a sin- 
gle constant scaling factor for the short-wavelength band and 
a factor with a linear dependence on frequency for the long- 
wavelength band. The resulting spectrum is shown in Fig. Q] 
(for clarity, we show the spectr um apodized using the extended 
Norton-Beer function 1.5 from N avlor & Ta hic 2007); we note 
that the short- and long-wavelength bands match very smoothly. 

Line fluxes were recovered from the calibrated unapodized 
spectrum using a custom-written tool. It first subtracts the 
underlying continuum using a grey-body fit, then it removes 
any remaining large-scale ripples using a polynomial func- 
tion. Emission lines were extracted by fitting sine-convolved 
Gaussian line profiles. The strongest line is fitted first and then 
subtracted, with the process repeated until no line is found above 
a pre-set discrimination level. The integrated line fluxes were ob- 
tained by calculating the area under the fit. Table Q] lists the line 
fluxes and their l<x uncertainties derived from the fitting proce- 
dure. In addition to the reported uncertainties we should include 
the following contributions: (i) the uncertainty in the estimation 
of the source-beam coupling factor due to the uncertainty in the 
beam profile and the assumption of identical distributions for 
dust and CO emission, and (ii) the uncertainty involved in the 
scaling the spectrum to match the photometric data, and in the 
measurement of photometric data. We conservatively suggest an 
uncertainty of ~30% for the line fluxes due to the above factors. 

Fig. E] shows the 12 CO SLED, which peaks at the J = 7-6 
line. In the same plot, we draw ground-based data compiled by 
Ward et al. (2003) (W03 hereafter). These data were measured 
with a smaller beam size, but were given for two observed lobes. 
In the plot we used the sum of the fluxes measured in the two 
lobes (which have a small overlap, but fit within the 43'.'4 beam). 



3. Radiative transfer modeling 



1 Extended-source flux calibration is derived from telescope emis- 
sion measurements, while the point-source flux calibration is based on 
observations of known astronomical point sources. 



We used RADEX (Ivan der Tak et al.l2007l) . a non-LTE code that 
computes the intensities of molecular lines by iteratively solving 
for statistical equilibrium using an escape-probability formalism 
assuming a uniform expanding sphere, to model the CO line in- 
tensities. The main inputs to RADEX are the gas density (n^), 
the kinetic temperature (7ki n ), and the CO column density per 
unit line width (A^( 12 CO)/Av). We ran the code for a large pa- 
rameter grid in Tidn (20-1 000 K), n H (10 2 -10 6 cm 4 ), JV( 12 CO) 
(10 I5 -10 18 cirT 2 ), andJV( 13 CO) (10 13 -10 I7 cirT 2 ). From this grid 
of models, we generated likelihood distributions by adapting the 
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Fig. 3. Left Panel: Likelihood distributions of kinetic temperature, density, CO column density and pressure. Right Panel: Likelihood 
contour plot of temperature and density. Dashed lines show constant pressure (Logio P (K cm -3 )) relations. 



Table 1. Measured fluxes of detected emission lines. Errors are 
lcr only from line fitting procedure, not including other uncer- 
tainties (see text). The 13 CO 7=6-5 line is missing due to fring- 
ing. 



Transition name 


Frequency 


Flux 






(rest, GHz) 


(10 3 Iykms- 1 ) 


(10- 16 Wnr 2 ) 


12 CO 7 = 4-3 


461.041 


74.1 ±2.2 


11.32 ±0.33 


12 CO 7 = 5^1 


576.268 


80.9 ± 2.3 


15.53 ± 0.45 


12 CO 7 = 6-5 


691.473 


74.0 ± 2.0 


17.04 ± 0.46 


l2 CO 7 = 7-6 


806.652 


77.7 ±3.1 


20.89 ± 0.84 


12 CO 7 = 8-7 


921.800 


60.7 ±2.1 


18.64 ± 0.65 


12 CO 7 = 9-8 


1036.912 


50.5 ± 2.3 


17.44 ± 0.79 


12 CO J = 10-9 


1151.985 


32.6 ± 1.3 


12.51 ±0.50 


12 coy = n-io 


1267.014 


21.9 ± 1.5 


9.28 ± 0.63 


12 coy= 12-n 


1381.995 


14.0 ± 1.2 


6.44 ± 0.57 


l2 COJ = 13-12 


1496.922 


7.1 ± 1.9 


3.53 ± 0.93 


13 CO 7 = 5-4 


550.926 


5.3 ± 0.7 


0.98 ±0.12 


13 CO 7 = 7-6 


771.184 


3.2 ± 0.6 


0.81 ±0.16 


"CO 7 = 8-7 


881.273 


2.3 ± 0.7 


0.68 ± 0.22 


HCN J = 6-5 


531.716 


2.9 ± 0.7 


0.52 ±0.12 


[Cl^P, -» 3 P 


492.161 


20.6 ± 1.6 


3.38 ± 0.26 


[CI] 3 P 2 ^ 3 P, 


809.342 


43.2 ± 0.9 


11.66 ±0.25 


[NII] 3 /> 1 -> 3 P 


1462.000 


124.1 ±5.8 


60.51 ±2.85 



method described in W03, for 7\ in , n Hl , 7Y( 12 CO), JV( 13 CO), and 
pressure by comparing the RADEX and observed line fluxes. 

To avoid any non-physical situation we applied two priors 
in this analysis following W03. The first one limits the 12 CO 
column density in a way that the total mass of the molecular gas 
producing the CO lines cannot exceed the dynamical mass of the 
galaxy according to the following relation: 



JV( 12 CO) < 



MdynXcO 1 



/j.m Hl 



< 



= 2.3x 10 20 cnT 2 



(1) 



where the dynamical mass of the disk Md yn = 2.4 x 10 9 M©, the 
radius of the disk R^ = 250 pc, the abundance of CO relative to 
H2, xco = 3 x 1CT 4 (W03), and ju = 1.4 is the mean molecular 
weight in units of wh, . The second prior limits the column length 
to be less than that of the entire molecular region according to 
^. (C : 0) < 1.54 x 10 21 cm. In this analysis we used all the 12 CO 

n(H 2 ).v C o J 

and 13 CO lines in TableQ]along with their lcr statistical errors. It 
was necessary to add 20% and 10% uncertainties for the CO 7 < 
8-7 and 7 > 8-7 lines, respectively, to avoid un-physically nar- 
row and noisy distributions (consistent with the additional 30% 
line flux uncertainty estimate in Sect. [2J. The resulting distribu- 
tions are shown in Fig. [3] for each variable, marginalizing over 
the other variables. The modeling only depends on the relative 
line fluxes, therefore the results will not be affected by the un- 
certainties in the absolute flux calibration. 
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Fig. 4. Comparing the highest likelihood model (dotted line) 
with our CO line intensities. The model shown by a dashed line 
was obtained by using only J > 7 -6 CO lines. The W03 data are 
shown by open circles. The inset highlights the deviations from 
models at the lower J end. 



4. Results and discussion 

We found that the highest likelihood model (dotted line in Fig. 
|4]i provides a good fit to our data (open squares), in particular for 
the higher 7 lines (7 > 7-6). The likelihood contour plot of tem- 
perature and density in Fig.[3](last panel) strongly indicates that 
the observed emission is coming from very warm gas with a ki- 
netic temperature of ~ 540 K and a pressure of ~ 3 x 10 6 K cm -3 . 
The detailed physical characteristics of the warm gas are listed 
in Table [2] which are obtained from the likelihood distributions 
shown in Fig. [3] 

The ISM of this galaxy has been well-studied using ground- 
based observations - in particular the lowest-lying CO rotational 
lines that provide constraints on the physical state of the cold 
molecular gas. Several studies from ground-based CO observa- 
tions, including W03, have identified cold gas at -30 K. We 
show W03 data (open circles) over-plotted in Figure [4] From 
7 < 7-6 lines, W03 deduced the presence of a warm component, 
in agreement with our finding, but with different temperature, 
mass and density. Having observations up to 7 = 13-12 enables 
a much better constraint to be placed on these parameters than is 
possible from the lower 7 lines. A colder gas component is also 
consistent with the deviation of our lower 7 lines, especially CO 
7 = 4-3, from our highest-likelihood model. If we use only the 
higher-7 lines (7 > 7-6) in our likelihood analysis, we get a 
model in better agreement with these lines while underpredict- 
ing the lower-7 lines, supporting the hypothesis of a contribution 
of the colder component to those lines. 

Assuming xqo of 3 x 10 -4 and an intrinsic line width of 180 
km s _1 (W03), and using our beam-averaged CO column density 
we find the mass of warm gas to be 1 .2 x 10 7 Mq within a beam 
area of about 2140 square arcseconds, likely covering most of 
the molecular emission from the galaxy. Using the mass of the 
cold gas from W03 we find a ratio of cold gas to the SPIRE- 
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Table 2. Model results and their uncertainties for the warm gas 



Quantity 



Most Probable Value range 1 



r kin (K) 


545 


350-825 


Log 10 n(H 2 ) (cm 3 ) 


3.7 


3.0-4.1 


Logi iV( 12 CO) (cm- 2 ) 


19.0 


18.5-19.8 


Log 10 © A Ar( 12 CO)*(cm- 2 ) 


17.4 


17.2-17.9 


7V( 12 CO)/7V( 13 CO) 


20 


15-37 


Log 10 P (K cm- 3 ) 


6.4 


5.8-6.7 


M gas (xl0 7 M o ) 


1.2 


0.7-3.6 



r Ranges are for 95% confidence intervals 

* Beam averaged column density where O a is an area filling factor 



observed warm gas mass of ~ 3. The best-fitting model predicts 
the optical depth for the CO lines, which peaks at a value of 1.7 
for / = 6-5, and then drops to approximately 10 2 for J — 13-12. 

Mid-IR H2 rotational lines are optically-thin and easily ther- 
malized, so they provide an independent constraint on the mass 
of warm gas. S everal transitions have been studied with ISO 
(Rigopoulouet al. 2002) and Spitzer (Beira o et al.l [2008). Both 
studies agree that the S(l) to S(2) line ratio suggests T ~ 500 K 
(assuming an ortho-to-para ratio of 3), in excellent agreement 
with our temperature. Using the Spitzer measurement of S(l) 
line flux corrected for our larger beam we calculate a mass of 
~ 2 x 1O 6 M0. Given the uncertainty on xco, and consider- 
able e xtinction (Td us t ^1-5; from dust models of lLaor & Draind 
(Il993l) extrapolated to our wavelengths) of the S(l) line, we find 
it to be consistent with our mass range. 

Our inferred thermal pressure (3 x 10 6 K crrr 3 ) is compara- 
ble to both that of the M82 ato mic gas as pr obed by the C II and 
OI transitions dKaufman et alJfl99 9: Lord et al.lll995l) . and the 
UV-shielded dense gas (Nayloret al. 2010), although this does 
not imply pressure equilibrium between the phases. Our warm- 
component mass is also similar to the 9 x 10 6 M0 inferred from 
the photo dissociation region (P DR) analysis based on the atomic 
gas lines (Kaufma n et al.lll999l) . However, the CO emission in 
the warm molecular gas likely does not arise from PDRs. This 
is because we measure CO transitions that are much more lumi- 
nous than predicted by the PDR models. These models require 
an extremely high density PDR (n > 10 5 ) to reproduce the J - 
6-5 and 7-6 intensities, a condition which is clearly ruled out by 
the atomic lines and their ratio to the far-IR flux. 

At a temperature of about 500 K, H2 will be the domi- 
nant coolant compared to CO. This is evident from the ob- 
served J^2jhTe_Juminosities, and agrees with the model predic- 
tions (lNeufeldetal.111995 : Le Bourlot et~ai1ll998l) . The models 



predict H2 cooling of ~ 10 ergs s per molecule for the 
temperature and density of SPIRE-observed warm gas. This im- 
plies a total molecular gas cooling of about 2.6 Lq/Mq; i.e., 
1.2 x 10 7 M Q of gas will radiate about 3 x 10 7 Lq in H 2 lines, 
in good agreement with the value derived from ISO and Spitzer 
data. 

What is the heating source of this warm molecular gas? 
Hard X-rays from an AG N have the po t ential for heating molec- 
ular gas in an XDR dMalonev et al.l 1 19961). but there is no 
strong evidence for an AGN in M82 (Strickland & Heckmanl 
2007). Moreove r, with a strong XDR component, such as 



seen in Mrk231 (Ivan der Werf et alJl2010h . the SLED becomes 
flat at high J instead of decreasing as in M82. Another sce- 
nario is heating via the enhanced cosmic ray density gen- 
erated by the high supernova rate in the nuclear starburst 
dSuchkov et al.l 119931) . With a cosmic ray energy deposition 
rate of 3.5-12 x 10 16 eV s -1 per H2 molecule in the Galaxy 



(iGoldsmith & Lan get! [1971 Ivan Dishoeck. & Blacld[l986l) . the 
energy input per mass in M82 is 0.09 to 0.3 Lq/Mq, too low to 
match the observed cooling. 

A second p ossibility i s heating from the dissipatio n of turbu- 
lence dFalgarone & Pugetll 19951: iMac LowHl999t iPan & Padoanl 
2009). Using a velocity gradient of 35 km s _1 pc -1 computed 
from our best fit RADEX model and a typical sizescale (or Jeans 
length) between 0.3 to 1 .6 pc in the expressio n for turbulent heat- 
ing per unit mass from Brad ford et al.l (2005), we can match the 
observed cooling of 2.6 Lq/Mq. Our velocity gradient is large, 
particularly when compared with the few km s 1 velocity spread 
found on 1 pc scales in Galactic star-forming sites, but may not 
be unreasonable in M82 given the powerful stellar winds known 
to be present in the starburst. 



5. Conclusions 

We have presented the Herschel-SPIRE spectroscopic observa- 
tions of the starburst galaxy M82. The spectra show a prominent 
CO emission-line ladder along with C I and N II lines. Radiative 
transfer modeling of CO lines clearly indicates a warm gas com- 
ponent at ~500 K in addition to the cold (~30 K) component 
derived from ground-based studies. The properties of the warm 
gas are strongly constrained by the high J lines, observed here 
for the first time. The temperature and mass of warm gas agree 
with the H2 rotational lines observations from Spitzer and ISO. 
At this temperature H2 is the dominant coolant instead of CO, 
and we argue that turbulence from stellar winds and supernovae 
may be the dominant heating mechanism. 
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